function tables_compstat_CCT(col_long,fp,output_dir)

% Comparative statics for CCT parameters:
% 1) baseline model
% 2) infinite CCT cost: kappa=infty
% 3) zero CCT cost but endo betas: kappa = 0
% 4) nonzero CCT cost but exogenous child discount factors: b2=b3=0
% 5) both zero CCT cost and exogenous betas: kappa=b2=b3=0

% Load appropriate results
load([output_dir 'temp_compstat_CCT'],'sim_data_all');

% Number of policy experiments (0 if we did no policies)
nr_exp = size(sim_data_all,4)-1;

% Number of models compared:
nr_models = size(sim_data_all,3); 
% Slice 1 = baseline

% Tables containing all the averages, for the baseline cases and the policy intervention cases
table_means = NaN(35,nr_models,nr_exp+1); % means across all ages

for i = 1:nr_models
     % loop over policy experiments
    for j = 1:(nr_exp+1)
        tmp_data = sim_data_all(:,:,i,j);
        
        % 5) CCT use, all ages 8-16
        tmprows = tmp_data(:,col_long.age) >= 8 & tmp_data(:,col_long.age) <= 16;
        tmp1 = tmp_data(tmprows, col_long.CCT); 
        table_means(5,i,j) = nanmean(tmp1);
        
        % 6) tau_c, all ages 8-16
        tmp1 = tmp_data(tmprows, col_long.tau_c); 
        table_means(6,i,j) = nanmean(tmp1);
        
        % 9) tau_p, age 8-16
        tmp1 = tmp_data(tmprows, col_long.tau1)+tmp_data(tmprows, col_long.tau2); 
        table_means(9,i,j) = nanmean(tmp1);
        
        % 23) Correlation between CCT and age_t, ages 8-16 (comparable to data moment)
        tmprows = tmp_data(:,col_long.age) >= 8 & tmp_data(:,col_long.age) <= 16;
        temp1 = tmp_data(tmprows,col_long.CCT);
        temp2 = tmp_data(tmprows,col_long.age);
        temp = nancov(temp1,temp2,1);
        table_means(23,i,j) = temp(1,2)...
            ./(nanstd(temp1(isnan(temp1) == 0 & isnan(temp2) == 0),1) ...
            .* nanstd(temp2(isnan(temp1) == 0 & isnan(temp2) == 0),1));
        
         % 24) Correlation between CCT and dad_educ, ages 8-16 (comparable to data moment)
        tmprows = tmp_data(:,col_long.age) >= 8 & tmp_data(:,col_long.age) <= 16;
        temp1 = tmp_data(tmprows,col_long.CCT);
        temp2 = tmp_data(tmprows,col_long.dad_educ);
        temp = nancov(temp1,temp2,1);
        table_means(24,i,j) = temp(1,2)...
            ./(nanstd(temp1(isnan(temp1) == 0 & isnan(temp2) == 0),1) ...
            .* nanstd(temp2(isnan(temp1) == 0 & isnan(temp2) == 0),1));
        
         % 25) Correlation between CCT and LW raw, ages 8-16 (comparable to data moment)
        tmprows = tmp_data(:,col_long.age) >= 8 & tmp_data(:,col_long.age) <= 16;
        temp1 = tmp_data(tmprows,col_long.CCT);
        temp2 = tmp_data(tmprows,col_long.lw_raw);
        temp = nancov(temp1,temp2,1);
        table_means(25,i,j) = temp(1,2)...
            ./(nanstd(temp1(isnan(temp1) == 0 & isnan(temp2) == 0),1) ...
            .* nanstd(temp2(isnan(temp1) == 0 & isnan(temp2) == 0),1));
       
        % 28) Final period child skills - ON GRID, log(k_17)
        tmp1 = tmp_data(tmp_data(:,col_long.age) == fp.M, col_long.log_k_17); % latent child capital at age fp.M+1 = 17 (ON grid values)
        table_means(28,i,j) = nanmean(tmp1);
        
        % 29) Final period child discount factor, betac_17
        tmp1 = tmp_data(tmp_data(:,col_long.age) == fp.M, col_long.betac_17); 
        table_means(29,i,j) = nanmean(tmp1);

        % 30) Final period child utility - u_{c,17}
        tmp1 = tmp_data(tmp_data(:,col_long.age) == fp.M, col_long.u_c_17); 
        table_means(30,i,j) = nanmean(tmp1);
        
        % 31) Final period child utility - u_{p,17}
        tmp1 = tmp_data(tmp_data(:,col_long.age) == fp.M, col_long.u_p_17); 
        table_means(31,i,j) = nanmean(tmp1);

        % 32-35) Average child/parents welfare W_c, w_p, discounted back to
        % the initial period (t_0, varies across hh), using the appropriate
        % discount factors for each household and age
        tmp_welfares = NaN(fp.R*fp.N,4); % col 1 = W_c,t0, col 2 = W_p,t0, col 3 = Bellman V_{c,t0} and col 4 = Bellman V_{p,t0}
        for r = 1:fp.R
            for q = 1:fp.N
                sim_id = r*1000 + q; % simulated id of this random household, see sim_choices_new.m to make sure this is the same
                beta_p = fp.beta_p_mat(q,r); % parents discount factor for this simulated household
                simdata_rq = tmp_data(tmp_data(:,col_long.id) == sim_id, :);
                seq_betac = simdata_rq(:,col_long.n)./(1+simdata_rq(:,col_long.n)); % 17 by 1 vector from ages 0 - 16
                seq_betac = [seq_betac; simdata_rq(fp.M+1,col_long.betac_17)]; % 18 by 1 vector, adding beta_c,17 at the bottom
                seq_uc = [simdata_rq(:,col_long.u_c); simdata_rq(fp.M+1,col_long.u_c_17) ];
                seq_up = [simdata_rq(:,col_long.u_p_tilde); simdata_rq(fp.M+1,col_long.u_p_17) ];
                
                tmp_t0 = simdata_rq(1,col_long.child_age97); % this household's initial age (so vectors start at index tmp_t0 + 1 and go up to 18)

                for t = fp.M:-1:tmp_t0
                    if t == fp.M
                        Wct = seq_uc(t+1,1) + seq_betac(t+1,1) * seq_uc(t+2,1); %W_16 = u_16 + beta_16 * u_17
                        Wpt = seq_up(t+1,1) + beta_p * seq_up(t+2,1); %W_16 = u_16 + beta_16 * u_17
                    else
                        Wct = seq_uc(t+1,1) + seq_betac(t+1,1) * Wct;
                        Wpt = seq_up(t+1,1) + beta_p * Wpt;
                    end
                end
                tmp_welfares((r-1)*fp.N+q,1) = Wct; % W_{c,t0}
                tmp_welfares((r-1)*fp.N+q,2) = Wpt; % W_{p,t0}
                tmp_welfares((r-1)*fp.N+q,3) = simdata_rq(tmp_t0+1,col_long.V_c);
                tmp_welfares((r-1)*fp.N+q,4) = simdata_rq(tmp_t0+1,col_long.V_p);
            end
        end
        
        table_means(32,i,j) = nanmean(tmp_welfares(:,1)); % W_{c,t0}
        table_means(33,i,j) = nanmean(tmp_welfares(:,2)); % W_{p,t0}
        table_means(34,i,j) = nanmean(tmp_welfares(:,3)); % V_{c,t0}
        table_means(35,i,j) = nanmean(tmp_welfares(:,4)); % V_{p,t0}
    end
end

% table_means

%% Now create a table with the summary statistics for the BASELINE models (no policy intervention)

tmp = [table_means(:,1,1) table_means(:,2,1) table_means(:,3,1) table_means(:,4,1) table_means(:,5,1)];

diary off;
diary_dir = [cd,'\table_compstat_CCT.txt'];
force_delete(diary_dir);
diary('table_compstat_CCT.txt');
diary on;

prec = '%.2f'; % precision for this table

disp(['\begin{table}[!htb]']);
disp(['\begin{center}']);
disp(['\caption{Comparative Statics - CCT Parameters}']);
disp(['\label{table:compstat_CCT}']);
disp(['\begin{tabular*}{\textwidth}{@{\extracolsep{\fill}}lccccc}']);
disp(['\hline \hline']);
disp([' & \multicolumn{1}{c}{\textbf{Baseline}} \rule{0pt}{3ex}',...
    '& \multicolumn{1}{c}{\textbf{(1)}} ',...
    '& \multicolumn{1}{c}{\textbf{(2)}} ',...
    '& \multicolumn{1}{c}{\textbf{(3)}} ',...
    '& \multicolumn{1}{c}{\textbf{(4)}} \\']);
disp(['\hline']);

disp([' \textbf{Model specifications} \rule{0pt}{3ex} & & & & \\']);
disp(['\quad - Util. Cost Mean/St.Dev. $\kappa$\rule{0pt}{3ex}          & 0.003 & $+\infty$    & 0         & 0.003    & 0 \\']);
disp(['\quad - Disc. factor process - $CCT$ slope $b^{up}_2$                  & -0.421 & -0.421       & -0.421    & 0         & 0 \\']);
disp(['\quad - Disc. factor process - $CCT \times age$ slope $b^{up}_3$       & -0.037 & -0.037       & -0.037    & 0         & 0 \\']);
disp(['\hline']);

disp([' \textbf{Simulated Moments} \rule{0pt}{3ex} & & & & \\']);

disp(['Fraction using CCT, ages 8-16 \rule{0pt}{3ex} & ',num2str(tmp(5,1),prec),' & ',num2str(tmp(5,2),prec),' & ',num2str(tmp(5,3),prec),' & ',num2str(tmp(5,4),prec),' & ',num2str(tmp(5,5),prec),' \\ ',...
    'Parental Investment Time $\tau_{p,t}$, ages 8-16 & ',num2str(tmp(9,1),prec),' & ',num2str(tmp(9,2),prec),' & ',num2str(tmp(9,3),prec),' & ',num2str(tmp(9,4),prec),' & ',num2str(tmp(9,5),prec),' \\ ',...
    'Self-investment Time $\tau_{c,t}$, ages 8-16 & ',num2str(tmp(6,1),prec),' & ',num2str(tmp(6,2),prec),' & ',num2str(tmp(6,3),prec),' & ',num2str(tmp(6,4),prec),' & ',num2str(tmp(6,5),prec),' \\ ',...
    'Corr. $(CCT_t,age_t)$, ages 8-16 \rule{0pt}{3ex} & ',num2str(tmp(23,1),prec),' & - & ',num2str(tmp(23,3),prec),' & ',num2str(tmp(23,4),prec),' & - \\ ',...
    'Corr. $(CCT_t,educ_2)$, ages 8-16 & ',num2str(tmp(24,1),prec),' & - & ',num2str(tmp(24,3),prec),' & ',num2str(tmp(24,4),prec),' & - \\ ',...
    'Corr. $(CCT_t,LW_t)$, ages 8-16 & ',num2str(tmp(25,1),prec),' & - & ',num2str(tmp(25,3),prec),' & ',num2str(tmp(25,4),prec),' & - \\ ']);

disp(['\textbf{Final-period Skills and Welfare} \rule{0pt}{3ex} & & & & \\']);

disp(['Final-period Child Quality $\log(k_{17})$ \rule{0pt}{3ex}& ',...
    num2str(tmp(28,1),prec),' & ',num2str(tmp(28,2),prec),' & ',num2str(tmp(28,3),prec),' & ',num2str(tmp(28,4),prec),' & ',num2str(tmp(28,5),prec),' \\ ', ...
    'Final-period Discount Factor $\beta_{c,17}$ & ',...
    num2str(tmp(29,1),prec),' & ',num2str(tmp(29,2),prec),' & ',num2str(tmp(29,3),prec),' & ',num2str(tmp(29,4),prec),' & ',num2str(tmp(29,5),prec),' \\ ']);
disp(['Final-period Child''s Utility $V_{c,17}(k_{17},\beta_{c,17})$ \rule{0pt}{3ex}& ',...
    num2str(tmp(30,1),prec),' & ',num2str(tmp(30,2),prec),' & ',num2str(tmp(30,3),prec),' & ',num2str(tmp(30,4),prec),' & ',num2str(tmp(30,5),prec),' \\ ', ...
    'Final-period Parents'' Utility $V_{p,17}(k_{17},\beta_{c,17})$ & ',...
    num2str(tmp(31,1),prec),' & ',num2str(tmp(31,2),prec),' & ',num2str(tmp(31,3),prec),' & ',num2str(tmp(31,4),prec),' & ',num2str(tmp(31,5),prec),' \\ ']);

disp(['\textbf{Other Welfare Outcomes} \rule{0pt}{3ex} & & & & \\']);

disp(['First-period Child''s Welfare $W_{c,t_0}$ \rule{0pt}{3ex}& ',...
    num2str(tmp(32,1),prec),' & ',num2str(tmp(32,2),prec),' & ',num2str(tmp(32,3),prec),' & ',num2str(tmp(32,4),prec),' & ',num2str(tmp(32,5),prec),' \\ ', ...
    'First-period Parents'' Welfare $W_{p,t_0}$ & ',...
    num2str(tmp(33,1),prec),' & ',num2str(tmp(33,2),prec),' & ',num2str(tmp(33,3),prec),' & ',num2str(tmp(33,4),prec),' & ',num2str(tmp(33,5),prec),' \\ ', ...
    'First-period Child''s Bellman Value $V_{c,t_0}$ \rule{0pt}{3ex}& ',...
    num2str(tmp(34,1),prec),' & ',num2str(tmp(34,2),prec),' & ',num2str(tmp(34,3),prec),' & ',num2str(tmp(34,4),prec),' & ',num2str(tmp(34,5),prec),' \\ ', ...
    'First-period Parents'' Bellman Value $V_{p,t_0}$ & ',...
    num2str(tmp(35,1),prec),' & ',num2str(tmp(35,2),prec),' & ',num2str(tmp(35,3),prec),' & ',num2str(tmp(35,4),prec),' & ',num2str(tmp(35,5),prec),' \\ ']);

disp(['\hline \hline']);
disp(['\end{tabular*}']);
disp(['\end{center}']);
disp(['\break']);
disp(['\textbf{Notes:} ``Baseline" corresponds to the baseline model. ',... %All values are averages across all households, time periods and simulation rounds. ',...    %'Final Child Quality ($\log(k_{17})$) is the latent value of child quality (in logs) at the end of age 16 or the start of period $t=17$. \\',...
    'Column (1) corresponds to the alternative model where CCTs are infinitely costly (i.e., $\kappa=\infty$). ',...
    'Column (2) considers the case where CCTs are costless (i.e., $\kappa=0$), but still endogenously affect the child''s future discount factor. ',...
    'Column (3) considers the case where CCTs are costly but no longer affect the child''s future discount factor (i.e., $b^{up}_2=b^{up}_3=0$). ',...
    'Column (4) considers the case where CCTs are costless and no longer affect the child''s future discount factor (i.e., $\kappa = b^{up}_2=b^{up}_3=0$). ',...
    'All experiments were done using $R = ',num2str(fp.R),'$ simulated data sets.']);
disp(['\end{table}']);

diary off;

end % function